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In this work we investigate the spectra of Laplacian matrices that determine many dynamic 
properties of scale-free networks below and at the percolation threshold. We use a replica formalism 
to develop analytically, based on an integral equation, a systematic way to determine the ensemble 
averaged eigenvalue spectrum for a general type of tree-like networks. Close to the percolation 
threshold we find characteristic scaling functions for the density of states p(A) of scale- free networks. 
p{\) shows characteristic power laws p(A) ~ A"^ or p(A) ~ A"^ for small A, where ai holds below 
and 02 at the percolation threshold. In the range where the spectra are accessible from a numerical 
diagonalization procedure the two methods lead to very similar results. 
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I. INTRODUCTION 

Recent studies of nets, ranging from social networks 
to power grids and the internet, revealed that in many 
cases the degree distribution pk, i.e. the probabil- 
ity that an arbitrary vertex is connected to exactly k 
other vertices, often exhibits a power law, namely that 
Pk ~ k~'' holds 1, 2]. Networks for which this rela- 
tion is fulfilled are called scale-free; scale-free networks 
differ from the classical random graphs for which 
the distribution pfcis Poissonian, and from small-world- 
networks 0) IE IE 0) 13 ■ Recent works have clarified 
that the properties of scale-free networks, in particu- 
lar percolation, differ markedly from the classical case 
Q II EE m El- It turns out that the asymptotic be- 
havior of Pk for k large, the so-called tail of pk, which 
is quantified by 7, is fundamental in differentiating be- 
tween the distinct classes of behavior: Thus for 7 < 4 the 
critical exponents change from the usual values found for 
classical graphs EjI • 

Now, the topological properties of a network are re- 
flected in the spectral properties of its connectivity ma- 
trix C: This matrix is constructed by letting its off- 
diagonal elements Cik be 1 if i and j are connected or 
otherwise; moreover, the diagonal elements Cu of C 
are zero. For scale-free networks it was found that the 
density of the eigenvalues of C has a triangular form 
with a power-law tail [13 . [TEl IT^ . On the other hand, 
many problems ranging from the dynamics of randomly 
branched polymers (IQj and the stress relaxation of near 
critical gels , over random resistor-capacitor networks 
I^H to glassy relaxation dynamics depend on the 
Laplacian A; A is connected to C via: 
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A whole series of works based on A were devoted to clas- 
sical deterministic and random graphs, such as Cayley- 
trees (dendrimers), hyper branched macromolecules, the 
Erdos-Renyi (ER) random graph and bond diluted 



Cayley-trees [Tl|Tl|Ti|2E|2l|2i|2l|2ll25','2^,'^|2i 
l2flj |. Based on the Laplacian, many time and frequency- 
dependent observables can be written as integrals over 
p{X), the density of eigenvalues of A: The structure of 
these observables is either 



Q{t)= / dXf{t,\)p{X), 



(2) 



or some related form, in which t is replaced by uj. For in- 
stance, for random walks, the site averaged return proba- 
bility Qnit) of a random walker to the origin is obtained 
with the choice /^(t. A) = e~'^* l2a|. Moreover, the me- 
chanical storage and loss moduli|2^ |2^ , the averaged 
time-dependent stretching of macromolecules in exter- 
nal fields [2E 113 , and the dielectric relaxation functions 
|30|. all obey forms similar to Eq. Q. We are inter- 
ested in the dynamic behavior of random graphs with 
arbitrary degree distributions and hence in the density 
p(A) of their eigenfrequencies. Following the ideas used 
in the analysis of gel dynamics poj ] and hyperbranched 
polymers [l9j, we display an integral equation for p(A) 
for a special class of random graphs with arbitrary de- 
gree distributions 0, |^ Is^, HE] ■ This integral equation 
allows us to determine p(X ) fo r the classes of scale-free 
networks discussed in Ref. [Tjl- 



II. RANDOM GRAPHS WITH ARBITRARY 
DEGREE DISTRIBUTIONS 

The ensemble of networks under consideration is ob- 
tained by starting from N vertices. Each vertex i has 
its degree fc^, and the probability distribution of the ki is 
Pk- As discussed in Ref. IIE, one can then connect the 
vertices pairwise through bonds (random pairing), while 
fulfilling the condition that the number of bonds ema- 
nating from each vertex i is given by its degree ki. All 
such possible combinations create the ensemble. In the 
limit N 00 the probability that a certain vertex is in- 
volved in a closed loop of bonds vanishes like 1/N [33l |: 
thus, in this limit a typical network realization is a set of 



2 



connected treelike clusters. Such a treelike structure may 
also be created as follows: We start from a vertex, say i, 
whose random degree k is chosen from the given distri- 
bution pk ■ Then each of the ki bonds of vertex i ends in 
a new vertex. One must note now that the probability of 
reaching via a randomly chosen bond a vertex of degree 
k is proportional to kpk, i.e. it obeys the distribution 



kpk 



EDO 
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It follows that we must now distribute the degrees of the 
newly produced vertices according to q^, Eq. (O. The 
procedure is then continued step after step and stops 
only where no new bonds were produced in the pre- 
vious step. The method creates all the random trees 
of the ensemble. Two examples of such ensembles are 
the bond diluted Cayley tree with functionality / [s^. 
where pk = {i)p'^{^ — p)'^~'^, and the Erdos-Renyi ran- 
dom graph, whose degree distribution is pk — p'^e^P /kl 



Depending on the p/j-distribution, the ensemble con- 
sists either exclusively of finite clusters or it includes an 
infinite, connected cluster, containing a finite fraction of 
the vertices of the system [s^. In Refs. 0, [s^l it was 
shown that the condition for the existence of this infinite 
cluster (also called percolating cluster) is given by 



> 1. 



(4) 



Eq. Q defines the so-called percolation threshold. It is 
very useful to extend the present model, by also allow- 
ing the strength of each bond to be weighted 19, 20] fol- 
lowing a normalized coupling strength distribution D{fi). 
Thus in the corresponding connectivity matrices, each of 
the nonzero values of Ctk can be chosen according to the 
distribution D{fi). 

As mentioned above, for a given network cluster S 
various dynamical quantities involve only the density of 
eigenvalues psW of the corresponding Laplacian A"^. 
Now, the ensemble averaged density of eigenvalues is 
given by 



pW = (ps(A)) =^wsps(A), 



(5) 



where the sum extends over all the clusters S, each of the 
/95(A) is normalized, and ws denotes the probability with 
which the cluster S is produced by the iterative growth 
procedure. Each of the S created in this way is con- 
nected, so that A'^ has only one zero eigenvalue, whose 



corresponding eigenvector is homogeneous. It turns out 
to be convenient to split off from p(A) the delta peak at 
A = whose weight is po, by setting: 



p{X)^p,5{X) + p+{X). 



(6) 



The density of states is connected with the diagonal el- 
ements of the resolvent, see e.g. Ref. Denoting 
for each site k the fc-th diagonal element of the resolvent 
R(A) = (A'5'-Al)"i byi?fcfc(A), what is needed is Rkk{X) 
averaged over all sites k and over all the ^-clusters: 



i?(A) = ((A^-Al)-i) 



(7) 



The probability Wk.s that when creating S we start at 
site k does not depend on fc; one has thus Wk^s = 'm^s/I'5'I, 
where \S\ denotes the number of vertices inside S. This 
leads to 



= T.T.^'^^siA' - ^^)kk 



S k=l 



(8) 
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Using for the normalized density of states of cluster S 
the relation 



PS I 



1^1 



fc=l 



we obtain from Eq. jSJ 



p(A) = lim -Imi?(A + ie), 



(10) 



with R{X) being given by Eq. I?)). Now, the average 
over the disorder can be performed with the help of the 
replica method [36| . 



III. DERIVATION OF THE INTEGRAL 
EQUATION 

In the following we denote the starting vertex by 0. To 
obtain the averaged trace of the resolvent we rewrite it 
with the help of a Gaussian integral over n-dimensional 
vectors 



RiX) 



Det- 
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see e.g. Ref. for details. The averaging procedure in Eq. Hll() is considerably simplified by taking the replica 
limit rt — > 0, since then the n/2-th power of the determinant is unity. Using 



(12) 



which follows readily from Eq. leads to 



i?(A)=- 
n 



i"o cxp 



4 s: '3 



cxp 



J<k 



(13) 



Here we use the dot over the equation sign to indicate 
that the limit n ^ has to be taken. Now we employ the 
fact that the S'-clusters are trees, in order to perform the 
integrations in Eq. p3fl iteratively, following the number 
g of growth steps. After the first growth step we have 



i?(i)(A)^- / dro r^exp 



n 

where we defined 
0(1) (ro) 

and 



.A 



Y^Pu{<i^^'Hv,)Y, (14) 



dr exp 
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i^(ro,r) 



F(rj,rfe)= / d^iD{^)cxp ~i^{rj-rk 



(15) 



(16) 



In a similar way, the averaged diagonal element after the 
second growth step reads 

/- ^ -1 OO 
dror^exp z^rg J2p>^{^^"^ (^o)}' (17) 
2 ^ — ' 
L J ; n 



with 



</.(2)(ro)= /drexp i^r^ F(ro, r) ^ g,{<^(i) (r)}'=-i. 

(18) 

More generally, introducing the generating functions of 
the pfe's and the g^'s: 



fe=0 k=l 

we find that after g growth steps 0(^-' (r) obeys: 
,A 



G^(l)' 
(19) 



(/)(») (ro) = y drexp 



F(ro,r)Gi[0(f-i)(r)], (20) 



and that it can be obtained iteratively, starting from 
0(1) (r) given by Eq. ifT^ . Furthermore, i?(f'(A) fulfills 



i?(»)(A)=- / dro r^exp 
n J 



•A 2 



Go[(/.(^)(ro)]. (21) 



Now, the n ^ limit can be performed as described in 
Ref. . This leads for 5 ^ oo to the pair of equations 



m = - 



A 



dxe'''Go[<j){x)] 



and 



(j>ix) = Oe--Gi[^{x)], 
where O is the linear operator 



(22) 



(23) 
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where (...)^ denotes the average over the distribution 

In Ref. [331 it was shown that the generating function 
Hoiz) = X^^i-fs^;* of the probabilities Pg that a ran- 
domly chosen vertex is part of a cluster of s vertices can 
be obtained based on the relations 

Ho{z) = zGo{Hi{z)) and Hi{z) = zGi{Hi{z)). (25) 

Here Hi{z) is the generating function for the distribution 
of sizes of components that are reached by choosing a 
random bond and following it to one of its ends. As a 
check for our scheme we now show that our Eqs. H22() 
and (|23|l are consistent with Eqs. H25|) . To do this we 
look for a solution of Eq. (|23|) in the form of a power 
series in A 



cj){x) 



E^' 

A;=0 



A>fe(a;). 



(26) 



We obtain the first term by comparing powers of A: 

0o(^) =e-"Gi(0o(a;)). (27) 

Identifying with z, Eq. H27() reproduces the second 
Eq. 12 5|) with 4)o{x) = Hi{e~^). Hence we infer from the 
first Eq. ^ that 



e-^Go(0o(x)) = Vp. 



s=l 



(28) 



4 



It follows that Go[(j)o{0)] is the probability for a vertex to 
be part of a finite size cluster. 

From Eqs. © and 1)10(1 one infers that i?(A) possesses 
a simple pole of the form po/^, where po is the finite 
weight of zero eigenvalues. Now po can be calculated by 
inserting Eq. (P?)]) into Eq. (|^ . which leads to 



Po = I dxe ^^GoiM^)) 



G'oil) 



(29) 



,(30) 



where in the second step we performed a partial integra- 
tion and used Eq. We note that inserting Eq. 
into (|30|l leads to po — Ps/s; the result represents 
the fact that each s-cluster contributes a term 1/s to the 
density of the eigenvalue zero. 



IV. SCALE-FREE NETWORKS CLOSE TO 
THEIR PERCOLATION THRESHOLD 

In this section we turn from our general considerations 
to focus on scale-free degree distributions; these exhibit 
for k large a power-law behavior, pk ^ k To describe 
the distance from the percolation threshold, Eq. we 
introduce the parameter A through the relation 



A = 1 



'EkLoHk - i)Pk 

T,T=0 ^Pk 



l-G'i(l). (31) 



Evidently, we assume by this that the first and the sec- 
ond moments of the pk distribution exist. From Eq. 
it follows that for A > the ensemble is made up of fi- 
nite connected clusters, while for A < there exists an 
infinite cluster. The critical point (percolation threshold) 
is at A = 0. As a note of caution we remark that the 
choice of the sign of A is possibly misleading but since in 
the following we investigate exclusively network clusters 
below the percolation threshold this choice considerably 
simplifies the formulas. Exemplarily, for the ER random 
graph the critical point is at Pc = 1 and for the bond 
diluted Cayley tree it is at = l/(/^ 1); using Eq. 
it turns out that in both cases A = [pc — p)/Pc- Note 
that for 7 < 3 Eq. 1(3 1() diverges; this agrees with the 
criterium of Eq. Q , since for 7 < 3 one always has an 
infinite cluster ilOj. 

We center now on the form of Eqs. ((23(1 and 1(22(1 close 
to A = 0. To be sufhciently general, we assume p^ to 
have for large k the form 



Pk 



k-''{c + 0{k''^)}. 



(32) 



This implies for 7 > 3 that the power series of Gi{4)), 
Eq. ifT^ . has as radius of convergence the unit circle 
|0| = 1. To determine the singularity on the radius of 
convergence we remark that the expansion coefficients 



Pk of the quantity Gi{(j)) = Gi{(j)) - cT{2 - 7)(1 - 0)^-^ 
obey pfc ~ k~^~^ for large k. Thus Gi(0) is m-times con- 
tinuously differentiable for |0| < 1, where m is the largest 
integer smaller than 7 — 2. Using the Taylor expansion 
of Gi (0) around t/i = 1 up to order m one gets ^2 113 



Gi(^) ~ 1 + (1-A)(,^-1) + 

+i^a;"Gi(i)(0-i)'^ 

777 ' ^ 



(33) 

-cr(2-7)(l- 0)^-2, 



where we used Eq. ((31(1 . Close to the percolation thresh- 
old A = we expect the solution of Eq. ((23(1 to scale 
in its variables x and A, and we choose a solution of the 
form [23 



(t){x) ~ 1 - A'30(a;/A^ A/A1+''), 



(34) 



with exponents /3 > and (5 > 0, to be determined below. 
Inserting Eq. ((34(1 into Eq. 1(231) and expanding in powers 
of A by using Eq. ((34(1 we obtain 



1- A-'x 



} 



l-A'30(x,A)= (35) 
[l-{p-')\^xdl + ...][ 

{1 - (1 - A)A'^0(x, A) + 0(7 - 4)^^A2'302(2;^ A) 
+^(4 - 7)cr(2 - -|)^^'^^'^^^-^{x, A) + . . . }, 

where 0{x) denotes the Heaviside function and the dots 
indicate terms with higher powers of A. Comparing pow- 
ers of A leads to the equation 



= {p-^)\/^^+l^xdl^{x) 
G'/(l) 



A^a 



A/3+1 



-^61(7-4)- 



+e{4: - 7)cr(2 - 7)A^('^-2)^7-2(2.)^ 



(36) 



Now the unknown exponents (3 and 5 are determined by 
the requirement that all terms in this equation be of the 
same order in A. For 7 > 4 this leads to (3 — 1 and 5 — 2. 
From Eq. it follows then: 



^ {^i-^)\xdl^{x) ~x + 4>{x) + 



(37) 



This universal scaling equation for the order parameter 
field (i>{x) was already pointed out in Ref. [2l| in con- 
nection with the mean-field theory of random resistor 
networks. Thus for 7 > 4 we obtain the classical mean- 
field scaling equation of the order parameter field 4'{x), 
which is also valid for the ER graph and the Cayley-tree. 
This is in accordance with the result in Ref. j!^], that 
the critical properties of classical random graphs are not 
changed for 7 > 4. On the other hand, for 3 < 7 < 4 the 
exponents (3 and S read now S ~ 1 + f3 and 5 ~ f3{j — 2), 
so that, solving for f3: 



3' 



(38) 
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as found in Refs. [Tl'| and Now the corresponding 
equation reads 

0= {^i-'^)Xxdl4>{x)-x + 4>{x) + cT{2--/)4>''~^{x). (39) 

We note that /3 is related to the probabihty Poo ~ A'^ 
that a vertex belongs to the percolating cluster. From 
Eqs. (|37|) and (|39|l for the order parameter field (p we 
obtain a scaling relation for R{X) by inserting Eq. I|34|) 
into Eq. (|22|l . To this end we subtract the pole po/X 
from -R(A) and expand in powers of A: 



R{X) 



Po 
X 



1 



A JO 

^l+/3 roo 



dxe-- {GoiH^)) - GoiMx))} 



X 



dx 



x{Go(l - A'^0(a:, A/A2+/5)) - Go(l - ^^M^))] 



A^- 



(fc) 



A/A 



2+13 



dx 



[4>{x,x/^''+p)~u^)]m 



where (^0(2:) = limA^o A"'^{(/io(a;Ai+'3) - 1}. Thus we 
have shown that p+(A) obeys for A ^ A close to a 
scaling law of the form 



p+(A,A)^A^-V~(A/A2+''). 



(41) 



Furthermore, the scaling function p{x) can be determined 
via Eqs. igUJl, 123, (EHl, and From the preceding 

considerations it follows immediately that the shape of 
p{X) differs in the region 7 > 4 from its shape in the 
region 3 < 7 < 4. In the first region Eq. H37|) is valid and 
p{X) does not depend on 7, whereas in the second region 
the 7-dependent Eq. holds. 



V. INTEGRATION FOR SCALE FREE DEGREE 
DISTRIBUTIONS AND SPECIAL 
DISTRIBUTIONS OF BOND STRENGTHS 

As shown in Ref. l20(, the analytical work simpli- 
fies considerably for the following distribution of bond 
strengths: 



D[p) = ^exp(-l/^), 



(42) 



since then the operator 6, Eq. H24I) . takes the form 



As noted in Ref. [23|, the particular choice of -D(^), Eq. 
(|42|1 . does not change much the small A behavior of p(A), 
given that in Eq. H42(l the probability for small coupling 
strengths p is exponentially small. In particular, D{p) 
docs not change the form of the function p(A), as only the 
first inverse moment (/i^^) enters Eqs. H37(l and (|39|l . Eq. 
(|44ll has to be solved subject to the boundary conditions 



(/)(0) 1 and (t){oo) = 0. 



(45) 



In the limit A ^ Eq. (|44|l can be linearized around the 
first term 4>[){x) of the asymptotic expansion, Eq. H26|l . 
This is achieved by inserting (j){x) = (f>o{x) + 4ii{x) into 
Eq. (|44|l and keeping only linear terms in 0;, since from 
Eq. if^ we have (/)i(a:) = 0(A). This results in the 
inhomogeneous linear equation 

Xxdlc^iix) + e-^G'Mx)]] Hx) = -Xxdl<t>o{x). 

(46) 

To investigate specific scale-free degree distributions 
we choose the following generating function 



A)(l 
1 



I7-3 



(1 



3(7-4)(7-2)(7-l) 



97-4 



(47) 



This form corresponds indeed to a degree distribution pk 
which obeys Eq. H32|) : the c- value in Eq. H32|l is 



(48) 



3(7-4)(7-2)(7-l)r(l-7)- 



The algebraic choice of Go(0) given by Eq. H47(l reduces 
considerably the effort needed to integrate Eq. H44|l . Fur- 
thermore, Eq. (|47() can be used in the whole interval 

3 < 7 < 5 containing the value 7 = 4 above which reg- 
ular mean-field exponents of percolation appear. On the 
other hand, for values of 7 outside the interval 3 < 7 < 5 
not all expansion coefficients of Go{(f>) are non- negative 
and thus they cannot be viewed anymore as probabili- 
ties. Note that the poles in 4 — 7 of the last two terms 
in Eq. (|47|l cancel and expanding Eq. (|47|l in powers of 

4 — 7 we obtain for 7 = 4 a branching point of the form 
(1 — 4>)^ log(l — 0) at = 1. Furthermore, from Eq. H30(l 
it follows that Go((/>) of Eq. leads to p+(A) being 
normalized to: 



dXp+{X) = i. 
^ 



(49) 



1 A 

6 = / d/x— exp(-l//i)exp xdl ^ [l + Xxdiy^ 

Jo M L J 

(43) 

For instance, applying 1 + Xxd^ = to both sides of 
Eq. H23|l . one obtains the ordinary second order differen- 
tial equation 



(x) + Xxd^^ix) ^ e-^'Gii^ix)]. 



(44) 



A. Numerical procedure 

To numerically calculate the eigenvalue spectra of 
scale-free networks we have performed extensive numeri- 
cal diagonalizations of the Laplacians of these structures. 

We create our structures by the recursive scheme in- 
troduced in the second section: For each realization of 
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the structure we first begin with an initial vertex, whose 
functionahty k is determined according to the probabil- 
ities pk derived from the generating function Eq. I|47|l . 
At the open end of each bond a new fc'-functional ver- 
tex is placed, now with k' distributed according to the 
probability distribution qk' given in Eq. iPJ; the latter 
procedure is then applied recursively to the k' — 1 open 
bonds of this new vertex. The recursion stops when no 
open bonds are left, i.e. when all outer vertices have func- 
tionality k' = 1. Note, however, that due to the limited 
time and memory resources available for the subsequent 
diagonalizations, the total number of bonds has to be 
restricted to some maximum value iVmax- If a given re- 
cursion has not stopped before reaching a total of NmnK 
bonds we proceed by closing all remaining open ends by 
a k' = 1-vertex and evaluate the properties of this trun- 
cated structure. Obviously, this also limits the range of 
validity of the resulting spectrum. As observed in our 
previous study in the region affected by the trun- 
cation the spectrum shows characteristic oscillations. To 
verify this procedure for the truncation limit A^max — 500 
used in general in this study, we have also performed for 
some of the curves shown in Fig. 4 additional diagonal- 
izations using A^max = 4000. The so obtained numerical 
results for logj^Q p{X) agree within the symbol size for the 
whole range covered by the A'max — 500 data shown in 
Fig. 4 and in fact extend the range of agreement with the 
solution of the differential Eq. (|44() by one order of mag- 
nitude. 

For the distribution of bond strength of a given struc- 
ture we chose either fixed bond strengths fi ~ I or 
strengths /i distributed according to Eq. (|42|l . The con- 
nectivity matrix Ay with entries weighted by these ^ is, 
by construction, a real, symmetric matrix. For all these 
matrices we obtained the eigenvalues using a combina- 
tion of the Householder method and of the tridiagonal 
QL diagonalization algorithm js^ [s^ • 

We accumulated the eigenvalues of all structures gen- 
erated for specific values of the parameters 7 and A, 
where eigenvalues stemming from a structure with \S\ 
monomers are weighted with a factor of 1/ 1 5*1, as given 
by Eqs. lO and For each of the data sets shown 
later in Figs 1, 3, and 4 the total number of structures 
truncated at A^max = 500 was 5 • lO*" and for structures 
truncated at A^max = 4000 was 4-10^. 




FIG. 1: Density of eigenvalues p(A) in double logarithmic 
scales. Displayed are spectra for the scale free degree distri- 
bution generated by Go{(t>), Eq. 1471 . and random coupling 
strengths /i obeying Eq. 11421 . Here, 7 = 3.5 is fixed and 
A is 0(0), 0.02(0), 0.04(V), 0.08(A), 0.16(o), and 0.32(0) 
from above. Lines: numerical solution of Eq. 1441 . Symbols: 
direct diagonalization of randomly created structures. 




-6 -5 -4 -3 -2 -1 
logioA 



FIG. 2: Density of eigenvalues p(A), obtained from the inte- 
gration of Eq. 1441 (right hand side of the figure) and from 
Eq. ignj (left hand side) for the same parameter values as in 
Fig0 To render this difi'erence clear, we have left out a small 
region around logj^Q A = —3.8. 



VI. RESULTS 

In Fig. n we display first the density of eigenvalues 
p(A) for the degree distribution pk generated by Go{(j)), 
Eq. I|47|l . with 7 = 3.5 and for various values of A < 0. 
The random coupling strengths p, obey the distribution 
of Eq. We obtained p(A) both through the numer- 

ical integration of Eq. (|44ll and also through the direct 
numerical diagonalization of many structure realizations, 
as described above. As can be seen, the p(A) obtained 
by the two methods agree very well with each other over 



a large A-range, thus supporting our theoretical consid- 
erations. The deviations of the curves from each other 
for small A are due to the limitations imposed by our 
direct diagonalization approach; in fact the sharp decay 
of the numerical results for A < 10"'^ is an artefact. The 
curves of Fig. ^ possess a shoulder at log^o A ~ —0.7 
which is most evident for the lowest curve corresponding 
to A = 0.32. However, this structure is caused by the 
choice Go(0), Eq. (|T7|) . and is not specific for scale- free 
degree distributions. As discussed above, scale-free net- 
works are characterized by the behavior for small values 
of A. To investigate the small A behavior, we show in Fig. 
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121 /o(A) for the same values of the parameters 7 and A, 
but extending to much smaller values of A. In this plot 
the results below A < lO^"^'^ are obtained from the lin- 
earized Eq. (|46|l . an approximation which we expect to 
be exact in the limit A 0. In fact, at A = 10~'^'^, where 
the approximation and the exact curve get together, the 
relative error of p(A) amounts to about 1%, which is al- 
ready hard to observe in the plots. For A not too close 
to the percolation threshold at A = and for small A 
the slopes of the curves in the double logarithmic plot of 
Fig. 121 tend to a constant. This would imply a simple 
algebraic dependence: 

p(A) -c(A)A"S for A^O, (50) 

with a positive exponent ai and a A-dependent coeffi- 
cient c(A). This differs from the situation expected to 
hold on classical random graphs with sufficiently fast de- 
caying degree distributions pk , where heuristic arguments 
have been given [2(1 123 | for the existence of Lifshitz tails 
in the density of states. In the latter situation one should 
observe the form nal 
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FIG. 3: Density of eigenvalues p(A) in double logarithmic 
scales at the percolation threshold A = for G(i{(f>), Eq. 14711 
and n obeying Eq. 14211 . Here 7 is varied, being taken to 
be 7 = 4.5(0), 4.25(0), 4(V), 3.75(a), 3.5(o) and 3.25(0) 
from above. Lines: integration of Eq. I|44|l . Symbols: direct 
diagonalization. 



p(A) - exp 



^(A) ' 

Va . 



for A ^ 0, 



(51) 



where A(A) ~ A'^/^ for A ^ 0. This behavior stems 
from the fact that small eigenvalues are produced by 
large, quasi linear regions, which however, occur with 
very small probability. Since for scale free degree distri- 
butions such linear regions are not likely to occur, there 
have to be other types of configurations which lead to 
an increase in the occurrence of small eigenvalues. For 
instance two vertices, each of very large degree, moving 
against each other produce a very low eigenvalue. 

At the percolation threshold A = we infer from Fig. 
2 for small A an algebraic decay of the form of Eq. H5()|l , 
with an exponent a2, which however differs from ai. One 
has namely a2 < o-i- Thus close to A = we encounter 
here a crossover behavior between two algebraic decays 
with different powers ai and a2- The scaling Eq. H41|) 
suggests that this crossover should take place at A ^ 

To determine the 7 dependence of a2, we display in 
Fig. inip(A) in double logarithmic scales for A = and for 
various values of 7. Assuming that for A — *■ the slopes 
of the plotted curves tend to a constant, say a2 we infer 
for 7 = 4.5, 4.25, 4, 3.75, 3.5 and 3.25 the values aa = 
0.017, 0.015, 0.04, 0.113, 0.25 and 0.5, respectively. For 
7 > 4 we expect to encounter the classical mean-field 
scaling function; since then p(A) tends to a constant for 
A ^ at A = [1^ l2^ , this implies that a2 = which is 
in good agreement with the first two values 02 — 0.017 for 
7 — 4.5 and a2 = 0.015 for 7 = 4.25. Moreover, directly 
at 7 = 4, we expect possible logarithmic corrections to 
the power law, Eq. H50|) which explains the value a2 = 
0.04 which is slightly too large. For 3 < 7 < 4 it turns out 
that we can reproduce the 7-dependence of 012 through 
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FIG. 4: Density of eigenvalues p(A) at the percolation thresh- 
old A = for fixed coupling strengths, /i = 1 in com- 
parison with the analytical results for fj, obeying Eq. H42|l 
(straight lines). The symbols show simulation data truncated 
at A^max = 500 (open symbols) and A^max = 4000 (filled sym- 
bols). The values of 7 and A and the symbol shapes are as 
in FigOl 



the relation: 



a2 



1 



27-5 



(52) 



In a similar way, we obtain from the small A behavior of 
p(A) for A > that ai = 27 - 5 holds. 

The 7 dependence of 02 can be derived from the scaling 
relation Eq. (|41|l by using the fact that the behavior of 
p(A) for A 3> 1 is given by the algebraic dependence of 



logio A 

FIG. 5: The quantity p(A, A) = A-V(AA*), for 7 = 3.5 in 
double logarithmic scales. The different curves correspond to 
A = 0.01, 0.02, 0.04, 0.08 and 0.16 from the right. 

p(A) for A < 1 at A = 

/5(A) = Uin A1-'^p(AA2+'3,A) - A1-'3(aa2+^)"^ 

For this to give a reasonable limit the exponent 1 ~ /3 + 
a2{2+/3) = of A has to vanish and solving this equation 
for a-z proves the relation Eq. H52I) . Similarly, the behav- 
ior of p(A) for A ^ 1 can be derived from the algebraic 
dependence Eq. ((SUJl of p{X) for A ^ 1 and A > 

p(A) = hm A1-''p(AA2+'5,A) - A1-'3c(A)(AA2+'5)"i 

= Ai-/3+"i(2+/3)c(A)A"i. (54) 

By the same arguments as above this leads to c(A) ~ 
A/3-i-ai(2+/3)^ but does not fix the exponent ai. 

As shown by Eqs. H39I) and H4U|) . the scaling function 
/5(A) of Eq. (|41|l should be model independent and there- 
fore should not depend on the particular choice of -D(/i), 
Eq. (|42|) . fn particular, at A = a power law with the 
exponent a2 should still hold. This is corroborated by 
Fig. ^ in which we display the density of eigenvalues 
p(A) for fixed bond strengths, ~ I, obtained from the 
direct diagonalization of random matrices, together with 
the analytical results for the distribution, Eq. I|42|l . of 
coupling strengths. 

To investigate the range of validity of the scaling law 
Eq. H41|) we display in double logarithmic scales for fixed 
7 = 3.5 in Fig. 5 and 7 = 4.5 in Fig. 6 the quantity 
/o(A,A) = Ai-^p(AA2+'3, A) which tends for A ^ to 
the scaling function /5(A) of Eq. H41|) . fn both figures 
we show curves for various values of A close to A = 0. 
For A and A small enough the curves for /5(A, A) should 
collapse into a single one, given by the scaling function 
/5(A) of Eq. (|41|l . In Fig. 5 the collapse appears roughly 
for AA-* < 10-2 and in Fig. 6 for AA^ < 10"^ 
Figs. 5 and 6 the scaling functions p(A) are given by 
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FIG. 6: The quantity p(A, A) = p(AA^) for 7 = 4.5 in double 
logarithmic scales. We show curves for A = 0.01 (diamonds), 
A = 0.02 (squares) and A — 0.04 (triangles). In addition, we 
display as a solid line the exact scaling function obtained by 
solving Eq. lIHTl . 



the envelopes of the curves /5(A, A) for different values of 
A. These envelopes p(A) seem to be monotonic growing 
functions. In Fig. 5 p(A) shows for the limiting case 
A the algebraic behavior /5(A) ~ A"^ but in Fig. 6 it 
is not possible to observe any algebraic dependence. 

In Fig. 6, the function p(A) given as a sohd hne is 
obtained from the direct integration of Eq. H37|) . which 
scales. More generally, for 7 > 4 we obtain f3 ~ 1; we 
note that for 7 > 4 all network ensembles lead to the 
same scaling function p{\). Note that the derivation of 
Eqs. (|37l) and H4U|) in Sec. IV shows that the essen- 
tial condition for this scaling function to hold is that the 
degree distribution pk decay faster than k~'^. This condi- 
tion is certainly fulfilled for classical random graphs, like 
the bond diluted Cayley tree or the Erdos-Renyi random 
graph. 



VII. CONCLUSIONS 

In this work we investigated the eigenvalues of Lapla- 
cians of structures belonging to a general type of tree-like 
networks, in which the vertex degrees are randomly dis- 
tributed. The Laplacian is of special interest, since it 
determines several, very important dynamic quantities 
associated with the network. For degree distributions 
Pk with a power law tail, pk ^ , we obtained p(A), 
the ensemble averaged density of eigenvalues, based on 
two different methods. First, in a traditional way, by 
performing numerical diagonalization techniques ,14 .ISj : 
second, using the replica method of statistical physics. 
The second approach allows to evaluate the ensemble av- 
eraged p(A) based on an analytical integral equation. For 
large A-domains it turns out that the agreement between 
the results obtained by the two methods is very good. 
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Of special interest is the behavior of p(A) close to the 
percolation threshold. Here an infinite cluster appears, 
and it is known that the exponent 7 which governs the 
large k behavior of pk affects the critical exponents of 
the percolation problem . With the help of our inte- 
gral equation approach we were able to study the scaling 
properties of p(A) close to the percolation threshold and 
to determine numerically the corresponding, 'y-dependent 
scaling functions. In agreement with Ref. |l2l |. we find 
that in the region 7 > 4 one recovers the critical proper- 
ties of classical random graphs. 

The long time dynamics is governed by the small A 
behavior of p(A). For this we found two algebraic forms 
p(A) ^ A"^ and p(A) ~ A"^ , where the first relation holds 
below and the second at the percolation threshold. On the 



basis of the numerical results of the integral equation we 
conjecture that ai = 2j — 5 and 0:2 = (4 — 7)/ai hold. 
We find that in scale-free networks very small eigenvalues 
occur with higher probability than in classical random 
graphs. We conjecture that this finding is due to the 
existence of highly connected vertices. 
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